#include "meshTypes.hpp"
#include "meshFunctions.hpp"
#include "meshTypeFunctions.hpp"
#include "vtkMeshExport.hpp"
#include "vtkMeshImport.hpp"
#include "util.hpp"
#include "streamAccess.hpp"

#include "testFunctions.hpp"
#include "mesh_testFunctions.hpp"

#include <catch2/catch_test_macros.hpp>
#include <catch2/catch_approx.hpp>
#include <sstream>

using Catch::Approx;
using namespace hms;

TEST_CASE( "Creating mesh from VTK XML" ){
	using namespace hms;
	using namespace hms::vtk;
	using namespace tinyxml2;
	MeshVariant meshVar;
	SECTION( "UniMesh" ){
		/* create mesh, then export it to vtk xml */
		meshVar.emplace<UniMesh>(
			4, 5, 0.1, 0.2, Vector2s{3.0, -0.2}
		);
		std::unique_ptr<xml::XMLDocument> meshDoc { meshToVtk(meshVar) };

		/* now build the mesh from that xml and check that they're the same */
		vtk::toMesh( meshDoc.get(), meshVar );

		REQUIRE( std::holds_alternative<UniMesh>(meshVar) );
		const UniMesh& mesh { std::get<UniMesh>(meshVar) };
		REQUIRE( mesh.nCellsX() == 4 );
		REQUIRE( mesh.nCellsY() == 5 );
		REQUIRE( mesh.dx() == Approx(0.1) );
		REQUIRE( mesh.dy() == Approx(0.2) );
		REQUIRE( mesh.origin().isApprox(Vector2s{3,-0.2}) );
	}
	SECTION( "RectMesh" ){
		ArrayXs xCoords (5);
		ArrayXs yCoords (7);
		xCoords << -0.1, 1.5, 3.0, 4.4, 5.7;
		yCoords <<  0.3, 0.5, 0.9, 1.5, 2.3, 3.3, 4.5;
		ArrayXs xCopy { xCoords.replicate(1,1) };
		ArrayXs yCopy { yCoords.replicate(1,1) };
		meshVar.emplace<RectMesh>(
			std::move(xCoords),
			std::move(yCoords)
		);
		std::unique_ptr<xml::XMLDocument> meshDoc { meshToVtk(meshVar) };

		vtk::toMesh( meshDoc.get(), meshVar );

		REQUIRE( std::holds_alternative<RectMesh>(meshVar) );
		const RectMesh& mesh { std::get<RectMesh>(meshVar) };
		REQUIRE( mesh.nCellsX() == 4 );
		REQUIRE( mesh.nCellsY() == 6 );
		REQUIRE( mesh.xCoordsInner().isApprox(xCopy) );
		REQUIRE( mesh.yCoordsInner().isApprox(yCopy) );
	}
	SECTION( "StructMesh" ){
		Matrix2Xs vertices (2,12);
		vertices <<
			0, 1, 3, 0, 1, 3, 0, 1, 3, 0, 1, 3,
			1, 1, 1, 3, 3, 3, 4, 4, 4, 8, 8, 8;
		/* distort the grid by translating every other vertex a bit */
		vertices(hms::all, Eigen::seqN(1,6,2)).colwise() += Vector2s{0.2,-0.3};
		Matrix2Xs vertsCopy { vertices(hms::all,hms::all) };
		meshVar.emplace<StructMesh>( std::move(vertices), 2, 3 );

		/* make some copies of data to be tested */
		const StructMesh& copy { std::get<StructMesh>(meshVar) };
		ArrayXs   lengthsCopy   { copy.edges().lengths().replicate(1,1) };
		Matrix2Xs normalsCopy   { copy.edges().normals().replicate(1,1) };
		Matrix2Xs centroidsCopy { copy.cells().centroids().replicate(1,1) };
		ArrayXs   areasCopy     { copy.cells().areas().replicate(1,1) };

		std::unique_ptr<xml::XMLDocument> meshDoc { meshToVtk(meshVar) };
		vtk::toMesh( meshDoc.get(), meshVar );

		REQUIRE( std::holds_alternative<StructMesh>(meshVar) );
		const StructMesh& mesh { std::get<StructMesh>(meshVar) };
		REQUIRE( mesh.vertices().isApprox(vertsCopy) );
		REQUIRE( mesh.nCellsX() == 2 );
		REQUIRE( mesh.nCellsY() == 3 );
		REQUIRE( mesh.edges().lengths().isApprox(lengthsCopy) );
		REQUIRE( mesh.edges().normals().isApprox(normalsCopy) );
		/* exclude corners which contain garbage */
		std::vector whichCells { 1,2,4,5,6,7,8,9,10,11,13,14 };
		REQUIRE( mesh.cells().centroids()(hms::all, whichCells).isApprox(
			centroidsCopy(hms::all, whichCells) )
		);
		REQUIRE( mesh.cells().areas()(whichCells).isApprox(
			areasCopy(whichCells) )
		);
	}

	SECTION( "Mesh" ){
		/* create triangle test mesh with the same algorithm for forming cells
		 * from vertices as is used in the importer */
		Mesh base { testMeshTriangles(true) };
		Array2Xs vertsCopy { base.vertices().array() };
		ArrayX<char> nDomainCellEdges { base.cells().nDomainCellEdges() };
		/* there are four cells with three edges each */
		Array4Xi domainCellVerts (4, 4);
		domainCellVerts = signalCode::uninitialised;
		for ( int i{0}; i<domainCellVerts.cols(); ++i ){
			domainCellVerts.col(i).head( nDomainCellEdges(i) ) =
				base.cells().vertices(i, true);
		}

		/* emplace test triangle mesh in variant */
		meshVar.emplace<Mesh>( cellVerticesToMesh(
			std::move(vertsCopy),
			std::move(nDomainCellEdges),
			std::move(domainCellVerts)
		) );

		/* now when exporting and reimporting the above mesh, all properties
		 * must be retained */
		std::unique_ptr<xml::XMLDocument> meshDoc { meshToVtk(meshVar) };
		MeshVariant reimported { Mesh{} };
		vtk::toMesh( meshDoc.get(), reimported );

		// std::cout << "Before:\n" << sstream.str() << "\nAfter:\n";
		// vtk::meshToVtk(std::cout, reimported);

		const Mesh& mesh { std::get<Mesh>(meshVar) };
		REQUIRE( std::holds_alternative<Mesh>(reimported) );
		const Mesh& copy { std::get<Mesh>(reimported) };


		REQUIRE( mesh.vertices().isApprox(copy.vertices()) );

		REQUIRE( mesh.nEdges() == copy.nEdges() );
		REQUIRE( mesh.nInnerEdges() == copy.nInnerEdges() );
		REQUIRE( mesh.nBoundaryEdges() == copy.nBoundaryEdges() );

		REQUIRE( mesh.nCells() == copy.nCells() );
		REQUIRE( mesh.nDomainCells() == copy.nDomainCells() );
		REQUIRE( mesh.nGhostCells() == copy.nGhostCells() );

		REQUIRE( mesh.edges().lengths().isApprox(copy.edges().lengths()) );
		REQUIRE( mesh.edges().normals().isApprox(copy.edges().normals()) );
		REQUIRE( ( mesh.edges().vertices() == copy.edges().vertices() ).all() );
		REQUIRE( ( mesh.edges().cells() == copy.edges().cells() ).all() );

		REQUIRE( mesh.cells().areas().isApprox(copy.cells().areas()) );
		REQUIRE( mesh.cells().centroids().isApprox(copy.cells().centroids(), hms::epsilon) );
		
		REQUIRE( ( mesh.cells().domainCellEdges().topRows<3>() == copy.cells().domainCellEdges().topRows<3>() ).all() );
		REQUIRE( ( mesh.cells().nDomainCellEdges() == copy.cells().nDomainCellEdges() ).all() );

		for (Index i{0}; i<mesh.nDomainCells(); ++i){
			REQUIRE( (mesh.cell(i).vertices() == copy.cell(i).vertices() ).all() );
		}
	}

	SECTION( "Mesh with boundary patches" ){
		meshVar.emplace<Mesh>( uniMeshAs<Mesh>(5,8,1.1,1.4,true) );
		MeshVariant importedMeshVar { Mesh{} };
		
		std::unique_ptr<xml::XMLDocument>
			meshDoc { meshToVtk(meshVar) },
			boundaryDoc { boundaryToXML(meshVar) };

		toMesh( meshDoc.get(), importedMeshVar );
		readBoundaryPatches( boundaryDoc.get(), importedMeshVar );


		REQUIRE( std::holds_alternative<Mesh>(importedMeshVar) );
		const Mesh& orig { std::get<Mesh>(meshVar) };
		const Mesh& copy { std::get<Mesh>(importedMeshVar) };

		REQUIRE( orig.vertices().isApprox(copy.vertices()) );

		REQUIRE( orig.nEdges() == copy.nEdges() );
		REQUIRE( orig.nInnerEdges() == copy.nInnerEdges() );
		REQUIRE( orig.nBoundaryEdges() == copy.nBoundaryEdges() );

		REQUIRE( orig.nCells() == copy.nCells() );
		REQUIRE( orig.nDomainCells() == copy.nDomainCells() );
		REQUIRE( orig.nGhostCells() == copy.nGhostCells() );
		
		REQUIRE( orig.edges().lengths().isApprox(copy.edges().lengths()) );
		REQUIRE( orig.edges().normals().isApprox(copy.edges().normals()) );
		REQUIRE( ( orig.edges().vertices() == copy.edges().vertices() ).all() );
		REQUIRE( ( orig.edges().cells() == copy.edges().cells() ).all() );

		REQUIRE( orig.cells().areas().isApprox(copy.cells().areas()) );
		REQUIRE( orig.cells().centroids().isApprox(copy.cells().centroids()) );
		REQUIRE( ( orig.cells().domainCellEdges() == copy.cells().domainCellEdges() ).all() );
		REQUIRE( ( orig.cells().nDomainCellEdges() == copy.cells().nDomainCellEdges() ).all() );

		for ( const BoundaryPatch& patch : orig.boundaryPatches() ){
			/* the patch with the same name must exist */
			REQUIRE_NOTHROW( copy.boundaryPatch( patch.name() ) );
			REQUIRE( (patch.edges()      == copy.boundaryPatch( patch.name() ).edges()      ).all() );
			REQUIRE( (patch.domainCells()== copy.boundaryPatch( patch.name() ).domainCells()).all() );
			REQUIRE( (patch.ghostCells() == copy.boundaryPatch( patch.name() ).ghostCells() ).all() );
		}

		// if (Mesh* mesh = std::get_if<Mesh>(&meshVar)){
		// 	std::cout << "Mesh with boundary patches:\n";
		// 	mesh->printInfo();
		// }
	}
}